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SYSTEM AVAILABILITY: TIME DEPENDENCE AND STATISTICAL INFERENCE 

BY (SEMI)NON-PARAMETRIC METHODS. 

D . P . Gave r 
P . A . J acobs 



INTRODUCTION. 

Component or system aval labi 1 ity generally reters to (a) the 
obability that the item oT concern is operable or ” up” and mission- 
pable at a point or during a period oT time, or to (b) the Traction 
time, or oT demands it the item is on standby, Tor which the system 
up. Applied to large systems such as entire U.S. Navy ships, system 
ailability is reTerred to as read i ness . At a diTTerent level, the 
ai lability and Tau 1 t-to 1 e ranee oT computer systems, both government 
d commercial, is oT concern, and has recently been extensively 
del led and simulated; see Goyal , Lavenberg and Trivedi . The 
ailability concept is also relevant when discussing the saTety and 
oductivity oT commercial nuclear power plants; in that arena it is 
antiTied by probabilistic risk assessment (PRA) . Related Tinite 
p.te stochastic models also occur in medical studies; c.T. Cox , 
cobs^ . 

Component or system availability is inTluenced both by the 
lerent Tai 1 ure-proneness oT the item and by the time and resources it 
<es to restore a Tailed item to service. Times to Tailure or ” up 
Ties”and to restoration or ’’down” times may vary considerably, and not 
cessarily independently, depending upon the mode oT Tailure, the time 
quired to diagnose the Tailure including the access to (including 
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competition i^or) diagnostic equipment and human skills, the 
availability of spare parts (logistic delays) , and other factors. It 
is also quite conceivable that the quality of the repair activity 
influences future times to failure. This effect is not recognized by 
the usual models; but see Gaver"^ (PP- 775-800), and Thomas, Jacobs and 
Gaver^ . Incorporation of ’’availability growth” in availability 
situations is a practical issue that is not frequently modelled; see 
our Section 2 below, however. 

The present paper addresses the assessment of simple system 
availability when there is concern about (a) t i me -dependence , so 
demands for system performance are not necessarily when the system is 
in ’’steady-state,” as is often assumed, and when (b) information about 
system failures and repairs is in the form of observed data so 
questions of statistical influence arise. The methods and models 
involved lean towards the sem i -paramet r i c or non-paramet r i c ; in 
particular we employ the empirical Laplace transform in the time- 
dependent scenarios of interest; non-paramet r i c estimation in assumed 
steady-state situations has been studied by Gaver and Chu® . Other 
investigations, e.g. Cox , Gross and Harris‘S and Ascher and Feingold^ 
have been overwhelmingly concerned with estimation in presumably well- 
specified stochastic (queueing, Markov) models; considerations of model 
m i s-spec if i cat i on are seldom broached. Here we propose analytically 
simple approximations to t i me -dependent system behavior, and assess the 
effects of model specification (”up” and ’’down” time dependence) upon 
rates of approach to a long-run steady state as the latter are 
estimated from available data (assumed to be a random sample). 
Uncertainty assessments (confidence limits and standard errors) are 
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furnished. More elaborate procedures involving Bayes or empirical 
Bayes setups that permit ’’strength borrowing” (in John Tukey^s phrase) 
are not addressed here, but are agenda items. 

□ur paper^s plan is as follows. Section 2 describes a selected 
group of probability models for simple system availability; no 
comprehensiveness is claimed. Solutions are given in terms of Laplace 
transforms, all of which are rendered immediately interpretable in 
terms of the random-t ime-of -demand or observation paradigm, described 
originally and applied to transform inversion in Gaver^®. Section 3 
describes and ’’fits” the simple exponent i al -approach -to-steady-state 
model used for representing t i me -dependent behavior; see Qdoni and 
Roth^ and Gaver and Jacobs where such an idea was used to represent 
time-dependent queueing behavior. Section 4 introduces issues of 
estimating t i me -dependent availability where only statistical data is 
at hand. Section 5 presents a variety of numerical illustrations. 
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2. TIME-DEPENDENT AVAILABILITY: MODELS, TRANSFORMS, AND THE RANDOM 

OBSERVATION (ROBS) INTERPRETATION. 

Consider a system that is active or operative (”up”) until Tailur 
occurs, after which a ’’down” period occurs, at the termination of whic 
a new up period begins, and so on. Suppose the system is to be 
responsive to a randomly appearing demand; what is the probability tha 
it will respond properly, or be available, at the time of the demand? 
In the military arena the device might be a surveillance (radar, sonar 
or communications system. In the commercial nuclear power context it 
might be a coolant pump. 

There are many versions of this problem. The simplest and most 
traditional imagines that the system operates continually between 
failures, and failure is instantly detected upon occurrence. Up times 
may be modelled as i id random variables, either directly or indirectly 
If the system is redundant a new up time begins at the instant an old 
one terminates, with repair off-line. Down times may also be iid. 
However, there are many plausible exceptions to these scenarios, a few 
of which are considered next. 

Here are various sample model formulations; the list is by no 
means complete. In all of these the system is up at t = 0; 
modifications are straightforward. 

Mode 1 1 : up times are iid {Uj} , with df Fy (x) ; down times are iid 

{Dj}, d.f. ^^dCy) • A('t) = P{System up at t | up time starting at t = 0}. 

Then, by a simple backward renewal argument, A(t) satisfies the 
integral equation 



A(t) 



Fu(t) 




A (t- v) F^ (d v) 



( 2 . 1 ) 
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where — F'(j*FQ(‘fc) , the convolution. Lapl a-ce-transForm i ng , 



A(s) = 



e’®^A(t)dt = 



1 - Fu(s) 

s[l - Fc(«)] 



( 2 . 2 ) 



If the item is observed or demanded at random time T ~ exp(s) then 



E[A(T)] = 



oo 






(2.3) 



and since E [T] = 1/s = r, t:he probability of being up upon demand is 



ao(T') 



E[A(T)] = |) 

1 - Fe(r-l) 



(2.4) 



E[U] 

E[C] 



E[U] 



E[U] + E[D] 



as r — 4. oo 



(2.5) 



by Taube r i an/Abe 1 i an results; see Feller (1966). One can 
the availability under random demand or observation; ROBS 
To emphasize the fact that exponential random observation 
we utilize EROBS . The last formula, written here as 



view (2.4) as 
for short . 
is involved 



a-o(oo) - E[U] + E[D] 

is the widely-applied long-run avai 1 ab i 1 i ty . 
priate in many of the contexts in which it is 
Note , though , that 



It may we 1 1 
app 1 i ed . 



( 2 . 6 ) 



be inappro- 






1 - E^(r-^) 

1 - Ec(^’') 



(2.7) 
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has a definite interpretation i^or al 1 r. It is usually easily compute 



and interpreted whenever transforms of U and D (or C) exist. As will 



be seen 5 it can be estimated from data under many c i rcumstances- -even 
when the model selected has either E [U] or E [D] infinite and (2.6) 
becomes uninteresting. 



directly written down, and directly interpreted under EROBS demands. 



i = 1,2,...} are iid and i nde pendent of Uq . Th i s attempts to mode 1 a 
situation in which the quality of maintenance during a down time 
affects the distribution of the next up time, so these times are 
dependent; the influence stops at that point, however. Let AQ(t) = 
P{system not up (down) at t | down time starts at t = 0}. Then 



Here are several additional models whose transforms can be 



Model 2 : 



The initial up time Uq has df Ey^ ; {Dj + Uj 








and 





AoCt-^Jo) = 1 - Ap(t-Uo) ^ t 



SO 




Transforming gives 
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1 - Fc(®) 



( 2 . 8 ) 



A(s) = i 



- ^UqCs) + Fuo(®)^d(s) 

1 - Fc(«) 



from which t-he observational probability aQ(r) can be immediately 
written down: aQ(r) = r‘^A(r’^). Models that propose that the up time 

can positively influence the duration of the following down time have 
been formulated; see Gaver^ . 

Mode 1 3 : (Random changepoint model; ’’reliability growth”). Suppose 

the system fails accord i ng to up t i mes {Uj^, i = 1,2,...!}, iid F-|^(u) 
and thereafter {Uj 2 > i = I + l , 1+2 , . . . } , iid p 2 (u); I ^ Pj i s the change- 
point and is associated with d i agnost i c- repai r activity. Let Fq(u) 
describe the iid down times {Dj , i = 1,2,...}. Suppose the system 
starts operation as a Type 1 up time, transitioning after a random 
number of down times to Type 2. Let Ajj(t) denote the probability that 
the system is up and in state j at time t given it is up at time 0, 
i.e. t is contained in an up time of Type j. Then, using backward 
renewal arguments , 



Aii(-fc) = Fiu(t) + (1-Pi) 



't 

Aii('b-r)Fi<-(dr) 

^0 



and 



A22(t) = F2u(t) + A22(t-r)F2c(di') , 

•'o 



Ai2(t) 



't 

A 22 (t-r)G(dr) , 

Jo 



(2.9) 



where 



G(r') = E PjF'ic(^) 5 

i=l 
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as before Fj^(r) refers fo fhe df of a cycle time Ujj + D j . Transform i 
immediately produces a solution: j 



All (s) 
A 22 (® ) 



1 ^ ~ ^iu(^) 

® 1 - (l-pi)Fic(s) 

1 ^ ~ ^2u(^) 

^ 1 “ ^2c(®) 



Ai2(®) — A22G(s) , 



(2.K 



and 

G(s) =_g Pj(Fic(s))‘ = p(Fic(s)) > 



this last being the generating function of the number of repairs to 
transition rv , I . Assemble to obtain the transform of the probability 
of being up at t: 




e (t ) 






Aii(s) + Ai2(s) 



1 

s 



1 - (1-Pi)^ic(®) 



+ p(^ic(s)) 



^ - ^2u(s) 

1 ~ F2c(®) 



( 2.11 



Consequently the observational probability 



ao(r) = E[A(T)] 



1 - Fnj(l/r) 

1 - (I-pOFicC^/"^) 



+ p(Fic(1/’’)) \ 



F2u(1/t) 

F2c(1/'^) 



(2.15 



E[U2] 

^ E[U2] + E[D] 



00 
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but the latter limiting formula provides no information about the 
availability at early times, possibly before change has taken place. 
The EROBS random demand formulation at least allows some information to 
be obtained in a very direct and meaningful manner by just evaluating 
the transform itself numerically. 

Mode 1 4 ; (Markovian changepo i nts ) . It is possible that a system 
alternates slowly between failure states, either temporarily or 
forever. The random motion may be in response to occasional changes in 
maintenance practices, to debugging (reliability growth), or to ageing 
(reliability decay) and subsequent replacement. To illustrate, let 
there be just two failure distribution states, as in Model 3 , but let 
{pjj ; i,j = be the transition probabilities of a Markov chain that 

governs jumps between them. Suppose the system starts in state 1 , at 
the beginning of an up time Uj ^ Fjy(r), i = 1 , 2 . Then, letting Ajj(t) 
be the probability the system is up and in state j at time t given it 
just started an up time at time 0 and was in state i. 



Aii(t) = Fiy(t) + pii [ Aii(t-r)Fi<-(dr) + f Aji (t-r) F^^- (dr) , 

J n J n 



•t ft 

= Pl 2 ^22 ("t-r) Fj^(- (dr) + p^^ Ai2("t-r )Fic(dr') , 

0 •'o 



( 2 . 13 ) 



ft ft 

l22(t) = F2u(t) + P22 A22(t-r)F2c(dr) + P21 Aj2 (t- r) F2C (d r) , 

J n J n 



't (*t 

A2i(‘t) = P21 Aii(t-r)F2c(dr) + P22 A2j (t- r) F2C (d r) 

Jq •’0 



Transforming yields "these equations, 
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Aii(s) — ^ ^ + Pll'^11 (®) ^IC (®) Pl2^21 (® ) ^IC (®) 

Ai2(s) = Pl2^^22 (®) ^IC (®) + Pll'^12 (®) ^IC (®) 

A22(s) = + P22A22(®)^2c(®) + P 2 I ^ 12 (® ) ^2C ( ® ) 

A2i(s) = P2iAii(s)F2c(s) + P22A2i(s)F2c(s) . 



Solve (2.14,a) and (2.14,d) simultaneously to get 



^11 



(s) = § 



^ ~ ^iu(^) 



1 - Pii^icC®) " Pl2P21^1C (®) ^2C (®) /(^ “P22^2C (®) ) 



From (2.14b,c), 



Ai 2 (s) 




1 - F^u(s) 

“P22^2C (®) ■"P2lPl2^1C (®) ^2C (®)/(^~Pll^lc(®)) 

/ Pl2^ic(^) ) 

^l-Pii^ic(s)^ 



It follows that if r = 1/s , 



ai(r) = sAit = s[Aii(s) + Ai2(s)] 



1 - Fiu(s) 

[1 “Pilule (®) “P 12 P 2 l^ic(®) ^ 2 C (®)/(^“P 22 ^ 2 c(®)) 

^ 1 ~ ^2u(^) 

1 -P22^2C (®)“Pl2P21^1c(®)^2c(®)/(^“Pll^lc(®)) 

( Pl2^ic(g) \ 



The above i^ormula is manii^est:ly opaque as it stands. 



Tau be r i an/ 



abelian results present the long-run behavior: 



ail (®o) 



1 i m sA 1 1 (s) 
s^O ^ 



E[Ui] 

E [Cl] + (P12/P21) E [C2] 



(2.18) 



tnd 



ai2 (00) 



E[U2] 

E[C2] + (P2i/Pi2)E[Ci] 



(2.19) 



>0 the sum 



Ll(oo) 



^11 (^) + ^12 (^) 



E[Ui] ^ E[U2] 

E[Ci] + (Pi2/P2i)E[C2] ^ E[C2] + (P2i/Pi2)E[CJ 

( 2 . 20 ) 



represents the long-run probability that the system is up. OT course 
L^(r) for r finite is available from (2.17) and potentially provides 
:onsiderably more information. 

Many other such models can be constructed and ’’solved” by 
transforms. In what follows we illustrate the way in which such 
transform solutions can be exploited to yield time-dependent and 
nferential information. The intention is to provide simple and 
'lexible approximate information rather than to utilize more exact, but 
lIso more computationally demanding, procedures. 
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3. TIME-DEPENDENT AVAILABILITY AT FIXED OBSERVATION TIMES (FOBS): 
SIMPLE EXPONENTIAL REPRESENTATIONS. 

II: is Familiar From Fhe Fheory oF Finite Markov chains and also 

From renewal theory that iF a long-run or steady-state condition is 
reached by a stochastic model, then the nature oF the approach is oFtei 
essentially exponential. In particular, iF, in Model 1, U ~ exp(A), 

D ^ exp(/i) then For A(0) = 1, 



A(t) = 



-f 



-(A-h/i)t 



A F /i A + /i 



(3.1) 



Such a time dependence is exhibited approximately by the M/G/1 queue; 
see Odoni and Roth^^; Morse^^’^^, and later Kielson^^, and others have 
assessed the relaxation times oF various queueing systems. See also 
Gaver and Jacobs^^ For a brieF discussion oF the time-dependent M/G/1 
queue utilizing such a representation. 

It is we I 1 -understood that natural generalizations oF the Markov 
Formulation producing (3.1) can. For certain models, lead to damped 
oscillatory approach to a steady-state value a(oo) , or to ultimate 
approach that is exponential modiFied by a negative power oF t. 
However, in what Follows we shall assume that the ’’error” at time t is 
approximately exponential: 



e(t) = A(t) - A(oo) 



ae 




(3. 2, a) 



or, more elaborately. 



(a^Fin2)o + (a^-ia2)e 



(3.2,b; 
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and i n ve s1: i ga1:e ways of assessing paramet:er values from (a.) perfectly 
specified models, wherein the distributions of U, D are presumed known, 
and (b) from data. This section deals with problem (a) . 

The motivation for considering simple exponential model (3.2) is 
the desire for an easily comprehended and computed assessment of time- 
dependent availability at a f i xed demand or observation time 
(abbreviated FOBS) . In many cases encountered such an assessment (or 
assessments) of oc and /? or of , Q 2 ? ^ and 9 can provide a useful sense 

of the behavior of A(t) as time progresses without the necessity of an 
extremely time-consuming and computer- i ntens i ve transform inversion 
procedure or of a symbolic inversion in terms of polynomial roots; the 
latter is ’’explicit,” but hardly comprehensible in general. 

Simulations and numerical solutions of governing integral equations are 
also useful approaches, but they are generally more computationally 
intensive than our proposals here. 

We suggest several ways of matching an exponential to A(t) by 
utilizing the transform A(s), presumed given. 

3.1 Method 1: Least Squares, Unweighted and Weighted (’’Tuned”) 

Begin by considering (3. 2, a). If 

e(t) = A(t) - A(c»); (3.3) 

we wish to represent this ’’error” by e (t) = ae ^ . Consider 

rOO rOO 

A(a,y0) = J ^e (t) -e (t) j^dt = J ^e ( t) -ae ^ j^dt , (3.4) 
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the integrated squared error; the object is to minimize A by choice o 
a and /? . Differentiation on a easily gives the optimizing condition 

Qq = [ e(t)e^^dt = 2/?e(/?) = 2(5 

i.e. is evaluated in terms of the known Laplace transform of A(t) 

and of e(t) , evaluated at (5, To find a minimizing (5 it is possible tc 
proceed by squaring the integrand of (3.4), discarding the e^(t) part, 
and substituting aQ for a. Simplification leads to 



A(/?) - 



A (oo) 

~T^ 






(3.5) 



-A(ao(/?) ,/?) = A(/?) = -2/?[e(/?)]2 



(3.6) 



which can be readily searched for a global maximum, /3q. It is easily 
verified analytically that if U ^ exp(A), D exp(/i) in Model 1 then 
Pq = X as it should. In general the search of (3.6) must be 

conducted numerically. 

The straightforward 1 east -squares procedure can be tuned towards 
different time ranges by appropriate weighting. We may seek to 
minimi ze 



A(a ,/? ; w) 




-ae ^^^^w(t;r)dt 



whe re 



w(t ; r) 



1 near r 
0 far from r. 



(3.7) 



(3.8) 
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len the minimizing a and f3 values, now orQ(r) , /?o(^) > should lead to an 

<ponential approximation that performs especially well near r — 
Ithough perhaps less well turther away. Again it is convenient to 
sleet an analytical torm tor w(t;r) that is compatible with the 
5 iplace transform, i.e. involves linear combinations ot exponentials, 
a example is the density tunction ot exponential order statistics; 
hi i s device was used in a numerical transform inversion scheme 
j?aver^^; Gaver and Jacobs^^ mention the idea in the present context) . 

Here is a specific example ot weighting. Let 




(3.9) 



here u = 7 ^. Then minimization leads T i rst to 
^ 6 r 




0 



hich gives 



e(/?+2/i) - e(/?+3/i) - 



2 /i 



1 




(3.10) 



r 



e(/?+2/i) - e(/?+3/i) 



~ l/(2/?+2/i) - l/(2/?+3/i) 



le((3+2fi) - e(/3+3/i)] 



(2/3+2//) (2/3+3/i) 



(3.11) 
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Squaring and omit:t:ing t:he term in e^(t) delivers 



-A(ao(/?) ,/?) = A(/?) = (2^+2/z) (2/?+3/x) [e(/?+2/i)-e(/?+3/z)]2 ^ 



which can be easily searched -for a global maximum. The estimation 



procedure is tuned towards r = in this case. A sharper tuning can 



be accomplished by a weight "function that concentrates more tightly 

around r than does (3.9); the latter is recognizably the density oT th 

median oT a sample oT size 3 "from an exponential d"f . IT w is the 

density oT the mth order statistic oT a sample oT n exponentials then 

its mean is 1 ^ ^ + ... -f ^ ^ ~ i ln(l - S) and its variance is 

n~l n-m-flj ^ 



Jk - > — — i-f m is proportional to n. As n increases the density 
^2 n(n-m+l) 

function oT w = w^ approaches a Gaussian/Normal density with variance 



decreasing like l/n, which is essentially a delta "function at 



T as n 



oo , 



3.2 Method 2: Derivative Matching: Exploiting Random Observation 
Times (ROBS) . 

Since it is proposed to approximate e(t) = A(t)-A(oo) by ae , an^ 
since the transforms o"f both sides are known (or can be estimated "fron 
data) one may contemplate determining a and /? as the solution o"f 



e(r) 



e‘‘^''^[A(u) -A(oo)] 4^ = iA (1/r) -A(oo) 

•*0 



pOO 



-u/r -0u 
e ae 



du 

T 



l+/^r 



( 3 . 13 ; 



that is, the system is observed at a random exponential time having 
mean r. One choice Tor r is t. Obviously any a, /? so found must 
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mplicitly depend upon r; in many situations for which the exponential 
s nearly correct the dependency of /? and a on r will be gentle. Here 
L,re two der i vat i ve- 1 ike prescriptions for finding a(r), /?(r): 

a) Solve 



5 + + A) = r(r + A) = l/e(r + A) 

5 + = l/e(^ - A) 



(3.14) 



or a = a(r;A), /? = /?(r;A). 

b) Invert and analytically differentiate at an interesting r-value 



r(r) 



T 

A ( 1 /r) - rA (oo) 




(3.15) 



r'(r) 



i 

a • 



(3.16) 



Now solve simultaneously for /?(r), a(r). A good diagnostic step is 

o plot Of(r) and /?(t) vs. r; if the plots are nearly constant with r 
hen the exponential provides a useful form; if the change is gradual 
t can potentially be represented by an empirical function, e.g. a low- 
rder polynomial . 

It is also possible to match the exponential approximation by r- 
uning, as was done by the weighting procedure of (3.7) . Simply 
bserve the system with probability density 



w(t;r) = c(n,m)(l-e ) e /i(e ) 



(3.17) 



t is convenient, to take n odd and m = [n/2] . 
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3.3 Met: hod 3: 



Ma1:ching Transform Means . 



The Tollowing procedure represent:s another possibility Tor 
assessing /? . Put a — 1 - A(oo) and write 



A(t) = A(oo) -f (1 - A(oo))P(t) 



(3.18 



and 5 denoting the transform of P by p> 




(3.19) 



se"^^ A (t ) dt 



sA (s) = A(oo) + (1 - A(oo))sp(s) 



Note that if 



A(t) 

1 - 



- A (oo) 
A (cxd) 




then transformation yields 



sp(s) 



1 

1 + nh 



(3.20 



o r 



[1 + /^(§)] sp(s) = 1 , 



(3.21 



SO 



_ 1 - sp(s) 

p(s) 
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To explicat:e the behavior of /? for large t , expand (3.22) in powers of 
s for s smal 1 : 



s [1 - sA (s) ] 

sA(s) - A(oo) 

sE[C]{Fc(s) - Fu(s)} 

E[C][1 - Fu(s)] - E[U][1 - Fc(s)] 



(3.23) 



sE[C]{-s(E[C]-E[U]) + ^(E[C^]-E[U^])+o(s^)> 
E[C]{sE[U] - ^E[U^]} - E[U]{sE[C] - ^E [C^] }+o (s^) 



Thus , 



lim /?(s) 
s^O 



2E[C]{E[C] - E[U]} 

E[C^]E[U] - E[C]E[U^] 



(3.24) 



It is anticipated that this expression is useful when t becomes large, 
provided a simple exponential is an appropriate approximation. 

We will evaluate (3.22) for s = 1 /t , t being the time of interest. 
It is also possible to tune this procedure to a time t as before. In 
particular, setting 




w(t ; r 




- A (oo) 
A (oo) 




pOO 

^0 



w (t ; r) e 



-/3t 



(3.25) 



and solving for (3 where, for example, w is as in (3.9) . This tactic 
will be employed for inferential purposes in Section 5 with w(t;r) of 
the form in (3.9) . 
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4 . TIME-DEPENDENT APPROACH CHARACTERISTICS ASSESSED FROM DATA. 

In This secTion we exa^mine siTuaTions rea.sonably represenTed as 
indicated previously, but Tor which limited inTormation on the up and 
down times is available: only samples oT Tinite sizes nu Tor the up 

times and n^ Tor the down times are at hand. Furthermore it is desirec 
to avoid using particular conventional analytical Torms Tor the 
underlying distributions; instead ” non-paramet r i c” or ’’distribution- 
Tree” methods will be invoked. Finally, the availability at a Tinite 
time aTter an initial moment is oT interest. 

4.1 Predicting Availability Assuming EROBS : Point and Interval 

Est i mates . 

Suppose we wish to estimate the availability oT a Model 1-type 
system at a random (exponential) time with mean r, i.e. assuming EROBS 
Data are available on the up times and down times: 

1 9 <^2 > • ♦ ♦ 9 

are the respective observations; although the order oT observation is 
Uj^ , d^ , U 2 5 d 2 5 . . . etc. , it will T i rst be assumed that the Uj^s and dj^s 
are independently sampled Trom Tixed distributions Fy and Fq , but the 
latter process is otherwise unspeciTied. 

A natural non - paramet r i c procedure is to deTine the sample or 
empirical transTorms 






nu 



nu 

E 



i=l 








( 4 . 1 ) 



and to utilize these in place of the population or true transforms thal 
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appear in (2.4) as a solution to the problem posed, 
parametric estimate of the desired availability is 



Thus a non- 



E[A(T)] = ao(r) = 



1 - 



1 - 






(4.2) 



Notice that the up and down times are conveniently assumed independent 
here. Procedures to deal with more complex models can and should be 
devised; accommodation to possible up and down data dependencies is a 
natural step, first by graphics and subsequently by model fitting and 
test i ng . 

There are various options for estimating the sampling variability 
of the estimate (4.2) given the basic model. A classical procedure, 
adopted here, is to first note that under our assumptions both <;^y(s) 
and <^q(s) are asymptotically Gaussian/Normal by the central limit 
theorem, for they are seen to be (modelled as) averages of independent 
and even bounded random variables e“^^ , e”^^ ; s = 1 /r > 0. We find 

easily that the corresponding random variables </>y(s), <^q(s) are 

approximately distributed as follows: 



. . ^ ^ Fu(2s) - (Fu(s))^ . 

~ N(^Fu(s), j ; 

, Fd(2s) - (Fd(s))2^ 



(4.3) 



Since a,Q(r) is a probability, confined between zero and one, it is 
reasonable to carry out further asymptotics on a transformation; we 
choose to study the logistic transformation 
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^(r) = In 



(^) 

l-ao(r) 



= In 



l-(^u(s) 



(4 



= In [l-«iy (s)] - In <^u(s) - ln[l-«iD(s)] 



Expand in Taylor series (t:he del1:a me1:hod) To obTain 



E[^(r)] In 



a-o(^) 



l-ao(r) 



a-o(^) 

l-ao(r) 



+ A 



1 Var['(u] + 1 - “g=^Pj 

^ (d-Fu)fuf ^ 



1 Var [<^p] 



+ i 



Sd/^d 



(4. 



where we have replaced Fy and Fq by Their esTimaTes <^y and (j>Q in The 
correcTion Term. FurThermore, 



Var[£(r)J ~ 



(Fu(I-^u))' 



Var[</>u] + -2 



®u / *^u 



_ „2 






2 “ 



(4 



We have put; 



I. 



and 



(4 






i=l 



Finally, approximaTe a-100% confidence limiTs are given by 
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C* + ■z.cti'i ■ Sg < In [ao(r)/(l-a.o(r))] < £* + z 



(4.8) 



3ind 






_ 



/ l-<Au(l/r) N 


if/ i-2<Au \ 


su . 1 ®dl 






(1-^d)" "-J 



(4.9) 



vh i 1 e S£ = from (4.6) . The expression (4.8) is Then inverted To 

provide Two-sided confidence limits for The acTuail ava. i 1 ab i 1 i Ty . A 
lumerical example is provided in a later section. 

4.2 Predicting Availability Sem i -Paramet r i cal 1 y at Fixed Time 
(According to FOBS). 

Utilize the same data as that in section 4.1 but suppose we wish 
bo estimate 



A(t) 2:^ A(oo) -f ae 



'/3t 



(4.10) 



\fter some experimentation it has been found that Method 3 above can be 
nost easily and effectively adopted to estimate /? and a: in (3.18) put 



A(oo) = _ 

u + d 



(4.11) 



i-nd estimate /? by evaluating the empirical version of (3.22). In the 
simplest case (Model 1) , the empirical version of p(s) is given by 



sp(s) = 



i-<;^u(s) _ u 

u + d 



(4.12) 



u + d 
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so the estimator ot (3.22) is 



3e(s) 



1 - sp(s) 

P(s) 



Sett i ng 



w(t;r) = 3/i(e - e 



in (3.25) results in the equation 



3/i[p(2/i) - p(3/i)] = 3/i 



(2// + /?) (3fi + /?) 



■] 



(4 . 13 



or 






H 

p(2/i) - p(3/i) 



0 



(4.14) 



a quadratic in /? . The positive solution to (4.14), gives another 

estimate oT /? : 



+ 4//[p(2/z) - p(3//)]‘^] (4.15) 

where /i = . In practice /3q may occasionally be negative or imagi- 

nary, in which case a reasonable alternative is required. In our simu 
lation, tests oT the procedure defaulted to an estimate appropriate fo 
U,D independent and exponentially distributed, i.e. ^ . How- 

ever in general a non-positive P may indicate inappropriateness of the 
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exponent, i al -approach model (4.10) so a problem-dependent alternative 
must be sought . 



The asymptotic behavior of 
here. In the following section 
use) procedures and results for 
estimates of A(t) and /? . 



the proposed estimates is not explored 
we describe bootstrapping (sample re- 
accessing sampling variability of 



f' 
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5. NUMERICAL ILLUSTRATIONS BY SIMULATION. 

The procedures described earlier will now be i 1 1 usL rated , and to 
some degree tested, using simulated data. 

5.1 Predicting Availability under EROBS : Illustrations 

ReTer to Section (4.1) Tor the basic approach to be utilized. 
Suppose that ny = n^ observations are available on i id exponentially 
distributed up and down time random variables; E [U] = 1, E [D] = 4. 

These data are to be processed utilizing expressions (4.8) and (4.9) . 
To do so, a total oT ns = 500 independent simulations were carried out 
Tor the cases ny = n^ = n = 10, and ny = n^ = n = 25. In Table 5.1 are 
reported the mean values oT the estimates oT availability at an 
exponential time with various means r. Results on 95% conTidence 
levels are also given (average values in parentheses) . 

Ill order that the point and interval estimates behave as well as 
they do Tor small samples (n « 10) the bias correction recognizable in 
(4.9) is required, and the normal percent-points are best replaced by 
Student t percent-points with n-1 degrees oT Treedom. An alternative 
to the above procedures is to apply the jackkniTe; see Gaver and Chu^ . 
However, the present method is perhaps more easily carried out on small 
computers utilizing nominal conTidence levels (l-a)‘100%. The Traction 
oT the conTidence limits covering/surrounding true availability at r 
was tabulated, as were the mean upper and lower conTidence limits on 
availability. The Tol lowing table provides a summary oT the results; 
as can be seen the coverage is close to the nominal 95%. 
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Tab 1 e 5.1 





Coverage and M 


ean Confidence 


Limits, EROBS 




Nom i nal 


Confidence Level 95% 
r- Val ues 




0.125 


1.00 


10.00 




True: (lower) 0.92 (upper) 


(lower) 0.83 (upper) 


(lower) 0.80 (upper) 


Availability: 


Estimated; (0.71) 0.90 (0.96) 


(0.66) 0.84 (0.93) 


(0.63) 0.81 (0.91) 


(n = 10) 


Coverage (%) 


94 


96 


94 


Availability: 


Estimated (0.82) 0.93 (0.97) 


(0.74) 0.83 (0.90) 


(0.70) 0.80 (0.88) 


(n = 25) 


Coverage (%) 


93 


97 


96 



The Tigures in parent:heses represent: mean upper and lower confidence 
limits; notice that, as anticipated, these tend to move towards the 
true values as the sample size increases. 

5 . 2 Predicting Availability under FOBS: Illustrations 

Again we use simulation to assess the accuracy of the proposed 
exponential approximations given in Section 3. 

It has been found experimentally that the emp i r i cal - t ransf orm- 
adapted least-squares approaches of Section 3.1 tend to be numerically 
unstable for small sample sizes; they behave well when component models 
are assumed perfectly specified (’’known’’) or for very large sample 
sizes. On the other hand, Method 3 of Section 3.3, adapted to the 
empirical transform as outlined in Section 4.2, performs satisfactorily 
for the various cases considered and is not computationally intensive. 

Simulation results for three different estimators for A(t) are 
reported. The estimators are all of the form 

A(t) = ^ -f (1 - ?)e (^*1) 
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where 



— i: =T 

u + d 

The estimators diTTer in the manner in which j3 is estimated. The T i rs 
estimate, Ap(t), estimates by 

/?P = (u-l + d-1) . (5.^ 

Note that Ap(t) is the estimator that would be obtained iT it were 
simply assumed that {Uj} and {Dj} are independent sequences oT 
independent exponential random variables and maximum likelihood is use 
to obtain estimates oT A(t) , i .e. iT the simplest Markov model were 
automatically invoked. The second estimator, Ae(t) , estimates (3 as in 
(4.13) with s = 1. The third estimator, Aq(t) , estimates as in 

pr 

(4.15) with /i = . IT the estimated ^ s are negative then P is set 

equal to /?p , the MLE estimator oT P iT it were known that {Uj} and {Dj} 
are independent sequences oT independent identically distributed 
exponential random variables. In our experience such pathological 
cases are rare . 

All simulations were done on an IBM 3033 AP at the Naval 

Postgraduate School using the LLRANDOMII random number generating 

17 

package; Lewis and Uribe^ . The simulations reported in Tables (5.2) 
and (5.3) have 500 replications. Generated in each replication is a 
sample of nu = 25 up times and n^ = 25 down times. In the experiment 
reported in Table (5.2) the up times and down times are independent. 
Three distributions are used to generate the down times: 
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Table 5.2 



Estimates of Availability 
Independent Up and Down Time.s 



Time Dist 


True 


/?P 




Ap 


(t) 


Pe 




Ae 


(t) 


/?q 




Aq(t) 






A(t) 


Mean 


Var 


Mean 


MSE 


Mean 


Var 


Mean 


MSE 


Mean 


Var 


Mean 


MSE 








(Median) 








(Median) 








(Median) 








.2 


A 


.83 


3.1 


.14 


.85 


.03 


3.8 


4.2 


.83 


.06 


3.8 


4.7 


.83 


.06 








(3.1) 








(3.4) 








(3.4) 








.5 


(Gamma) 


.71 


3.1 


.14 


.74 


.05 


4.6 


34 


.72 


.06 


4.5 


23 


.72 


.06 








(3.1) 








(3.6) 








(3.7) 








1.0 




.67 


3.1 


.14 


.68 


.05 


4.7 


48 


.68 


.05 


4.6 


23 


.68 


.05 








(3.1) 








(3.6) 








(3.5) 








.2 


B 


.85 


3.1 


.23 


.85 


.03 


3.4 


3.2 


.85 


.05 


3.3 


3.1 


.85 


.06 








(3.1) 








(3.0) 








(3.0) 








.5 


(Exp) 


.74 


3.1 


.23 


.74 


.04 


3.5 


4.7 


.74 


.06 


3.5 


4.7 


.74 


.06 








(3.1) 








(3.0) 








(3.0) 








1.0 




.68 


3.1 


.23 


.68 


.05 


3.7 


14 


.69 


.06 


3.9 


18 


.69 


.06 








(3.1) 








(3.0) 








(3.0) 









.2 


C .90 


4.4 

(3.4) 


6.7 


.86 


.06 


3.4 

(1.7) 


38 


.91 


.03 


3.5 

(1.8) 


52 


.91 


.03 


.5 


(M. Exp.) .87 


4.4 

(3.4) 


6.7 


.77 


.13 


2.9 

(1.2) 


57 


.87 


.04 


2.5 

(1.2) 


14 


.87 


.04 


1.0 


.84 


4.4 

(3.4) 


6.7 


.72 


.17 


2.5 

(1.0) 


26 


.82 


.06 


1.9 

(.92) 


6.1 


.83 


.06 


2.0 


.79 


4.4 

(3.4) 


6.7 


.70 


.17 


2.5 

(.87) 


75 


.77 


.09 


1.8 

(.80) 


9.8 


.78 


.08 











Table 5 . 3 
















Est images 


of Availability 














Dependent: 


Up and Down 


T i mes 




















•Ac 




(b) 


T i me 


Case 


T rue 






Aq ( t ) 


/?q 




Aq 






A(t) 


Mean 


Var 


Mean MSE 


Mean 


Var 


Mean 


MSE 


0.2 


DA 


.94 


3.1 


1 .7 


.85 .10 


1 . 0 


. 13 


.94 


.02 


0.5 




.87 


3.1 


1 . 1 


.75 .13 


1 .0 


.07 


.87 


.03 


1 .0 




.79 


3.1 


1 .0 


.69 .11 


1 . 0 


.07 


.79 


.03 


2.0 




.72 


3.2 


1 .9 


.67 .05 


1 . 1 


.09 


.71 


.02 


0.2 


DB 


.97 


3.1 


14 


.91 .07 


.73 


5 . 1 


.97 


.02 


0.5 




.93 


2.4 


11 


.87 .07 


.83 


3.8 


.93 


.03 


1 .0 




.88 


2.0 


7.2 


.83 .07 


.79 


.56 


.88 


. 05 


2.0 




.81 


1 .8 


8.0 


.78 .08 


.83 


.45 


.81 


.07 
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A: 


P{D > -t} 


= 1 - - dte'"*^ , t > 


0 ; (Gamma) 




B: 


a 

V 

ct 


= e'^‘ , t > 0 ; 


(Exponent i al ) 




C: 


P{D > -t} 


= .9e‘®^ + t > 0 ; 


(Mixed Exponential) . 




The 


up times are generated Trom an exponential distribution with 


unit 


mean 


. The theoretical values oT A(t) 


in each oT these cases are 




A: 


A(t) = . 


67 + e''*-®‘[.33 cos(1.9t.) 


-f . 26 si n (1 . 9t)] ; 


(5.3) 


B: 


A(t) = . 


67 + .33e'^‘ ; 




(5.4) 


C: 


A(t) = . 


67 + .24e’-^'‘‘ + .09e‘®'^‘ . 




(5.5) 



For each replication ot the simulations whose results are reported in 
Table (5.2) the three estimates oT A(t) are computed. The mean, and 
mean square error (MSE) , oT Ap(t), Ae(t) and Aq(t) are computed: i.e. 

the mean is 

1 500 . 

= ^k-i 

and 

1 500 _ - 

" A(t))2 (5.7) 

where Ae(k;t) is the point estimate at t in the kth realization. 
Furthermore, the mean and variance oT the estimates oT P are displayed 

5 . 3 Discussion oT Tables . 

Table (5.2) reports results Tor simulations with data sample size 
25. Not surprisingly, the exponential estimator Ap(t) has means equal 
to the values of A(t) and the smallest mean square error in the case B 
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exponent:ia-l down *times. The meains oT Ap(t:) are somewhat: diTferent. 
Trom the A(t)^s in case A (gamma down times) and quite ditferent in 
case C (mixed exponential down times). The estimators Ae(t) and Aq(t) 
have means which are closer to the A(t)^s than those of Ap(t) in cases 
A and C. The means of Ae(t) and Aq(t) are a bit low for t = 1, 2 in 

case C. In cases A and B the values of the MSE of Ae(t) and Aq(t) 
indicate the greater variability of these estimators. In the 
exponential case B, the true value of /? is 3; the estimator ^p has 
means closest to 3 for case B and much smaller variances compared to 
and /?q . In the other two cases p(t) is not of the form e ; the 

theoretical values of A(t) are given in (5.3) and (5.5). 

In case A (gamma distributed down times) the single parameter, /? , 
of exponential decay is 4.5. The means of the estimators, and /?q are 

closer to this value than those of /?p . In case C of mixed exponential 
down times, the smallest parameter of exponential decay is 0.341. Once 
again the means of a.nd /?q are closer to this value than are those of 
/3p . The variances of /?q and /?g for the sample size of 25 can be large. 
In many cases the variance of /?q is less than that for /?g • Increasing 
the sample size decreases the variances, as is anticipated. 

Table (5.3) summarizes a simulation study of the procedure of 
Section 4.2 in two cases in which the pairs (Uj,Dj) are independently 
and identically distributed with Uj and Dj perfectly dependent. The two 
cases are : 

DA: Uj = Ej and D| = ^Ej with {Ej} independent identically 

distributed exponentials with unit mean. 

{ ^Ej with probability .9, 

4E| with probability .1 



31 



with {Ej} independent exponential random variables with unit mean, 
case DA 



A("t) = § + 5 • (5.8 

In case DB 

A(t) = § + 0.1 le'^ + 0.22e'°-^^* . (5.9; 



In the simulations two estimators tor F^(s) are used. The First 
(incorrectly) assumes that {Uj} and {Dj} are independent and uses 
(^u(s)<^d(®) estimate F^(s) . The second does not assume independence, 

and instead estimates F^(s) by 






nr ^ ^ 

^ i=l 



(5.10 



The table shows results tor the estimators Aq(t) and ySq . The sample 
size is 25 and there are 500 replications. 

Table (5.3) intorms us that the estimators may be noticeably 
sensitive to the choice ot the estimator ot F^(s) . It Uj and D| are 
dependent, then using <^u(®)^d(®) estimate F^(s) can be quite 

misleading. As a result, it is suggested that it there is a 
possibility that Uj and Dj may be dependent, 0^(s) is the more model- 
robust estimate ot F^(s), and hence ot the desired availabilities. 

Note that in case DA, A(t) has exactly the presumed torm ot 
(4.10). When 0^ is used in this case, the mean value ot the ySq ^ s equal 
tlie theoretical value ot unity in all but the case t = 2.0. In case 
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DB 5 p(t:) is actually a mixture ot two exponentials, one having rate 

unity and the other having rate 0.27; when <f>Q is being used, the mean 
values of the ^q^s fall between the two correct rates. Emphasis of 
either rate, and the corresponding probability, can be achieved by 
we ight i ng . 

5 . 4 Sem i -Paramet r i c Confidence Intervals Illustrated . 

Table (5.4) shows bootstrap confidence intervals for Aq(t) and /?q 
at various times. A single sample of 25 up times and 25 down times is 
generated by simulation. The up times are independent with unit-mean 
exponential distribution; they are independent of the down times. The 
down times are independent with the gamma distribution, case B. One 
hundred bootstrap replications were then carried out and Aq and are 

computed for each replication. Both estimators of F^(s), (s) <^q (s) 

and <;6^(s), are used. Table (5.4) exhibits the 5th, 50th, and 95th 
order statistics of the estimates which give 95% confidence intervals 
of the parameters. 

All confidence intervals for /?q cover the true values for the 
particular sample utilized. The confidence intervals for /?q cover the 
exponential decay parameter of value 4.5, for the gamma case. The 
confidence intervals for Aq(t) are the same for both methods of 
estimating F^(s) suggesting that not much is lost by using <?^^(s) to 
estimate F^(s) in the independent case. As a result it is suggested 
that if there is a possibility that {Uj and {Dj} may be dependent then 
<^^(s) is the more robust estimate of F^(s) . 

A similar bootstrap experiment was carried out with a single 
sample of 25 up times and 25 down times generated from model DA in 
which the up and down times are dependent. One hundred bootstrap 
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Table 5.4 



Bootstrap ConTidence Intervals 
Gamma Distributed Down Times 
Independent Up and Down Times 



























Time 


True 






q 




Aq(t) 






q 




Aql 


(t) 




A(t) 


.05 


.50 


.95 


.05 


.50 


.95 


.05 


.50 


.95 


.05 


.50 


.95 


.2 


.83 


1.6 


4.1 


13.6 


.71 


.84 


.92 


1.6 


3.9 


11.8 


.71 


.84 


.92 


.5 


.71 


2.4 


4.1 


12.9 


.61 


.73 


.82 


2.4 


3.9 


9.8 


.61 


.73 


.82 


1 


.67 


2.3 


3.4 


16.3 


.58 


.71 


.79 


2.2 


3.4 


19.1 


.58 


.71 


.79 



replications were generated and Aq and computed Tor each replication 

with both estimators oT F^~(s). Once again the resulting conTidence 
intervals Tor A(t) were very similar Tor both estimators oT F(-(s) . 

OT course the above results are quite Tragmentary, but seem useTul 
and promising. Further sampling experiments and asymptotic analyses 
will shed more light on the behavior oT the estimating procedures 
explored, and may well suggest alternations or replacements. 
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6. CONCLUDING COMMENTS. 



This paper argues that: probability models oT the availability of 
various systems can be expressed in terms of Laplace transforms, and 
that the Tinite-time behavior of such systems can be inferred, nearly 
non-paramet r i cal 1 y , Trom data. Our approach has been to invoke the 
empirical Laplace transforms and to utilize its easy direct interpreta- 
tion (EROBS) , in conjunction with a presumed approximate exponential 
rate of approach to steady state, to deduce availability at a fixed 
finite time after a known initial moment (FOBS). 

The methods proposed are distinguished by their simplicity and 
moderate computer intensivity as well as by their lack of direct 
dependence upon probability models in “up” and “down” times selected 
from conventional families such as the Gamma. Simulations have been 
used to evaluate the procedures suggested, and to provide approximate 
confidence limits, either by asumptotics (utilizing the approximately 
Normal /Gauss i an behavior of the empirical transform), or by a simple 
re-sampling, Efron ^s bootstrap . In particular, we have examined the 
effect on inference quality of assuming the wrong j o i nt probability 
model: one that falaciously assumes independence when dependence 

(between up and down times) is actually present. 

No claim is made that the methods proposed here are the best 
available; in fact there are many alternatives. One is to analytically 
invert the empirical transform of availability possibly by use of the 
Stehfert algorithm^^, although competitors are available. Another is 
to bootstrap directly; the latter exercise involves re-sampl ing up and 
down times from the observed data, to reconstruct the sample path of 
the process, and to score 1 at time t if an up time covers t, otherwise 
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score 0. Finally, A(t:) is estimated by the proportion of (re)samples 
that count 1. Confidence limits are available from the basic bootstra 
technology. The direct bootstrap approach is being investigated by 
Lee'^^ . Another option is to replace the distribution functions in the 
renewal equations for availability by their empirical counterparts and 
numerically solve the empirical renewal ( Vo 1 te r ra-type ) equations, wit 
subsequent bootstrap follow-ups to assess uncertainty. All such 
methods promise to be far more computationally intensive than our 
present approximate approaches. Their investigation has been deferred 
Application of our approximation procedure to infer the M/G/1 queue 



finite-time behavior is under way; Jacobs and Gaver 
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